Dataset ejemplo: lung (paquete
survival)
#| echo: true
library(survival)
data(cancer, package="survival")
table(lung$ph.ecog)
0 1 2 3
63 113 50 1
lung <- lung %>%
na.omit() %>%
mutate(status2 = ifelse(status == 2, 1, 0)) %>%
filter(ph.ecog<3) %>%
mutate(surv_obj = Surv(time, status2),
sex = factor(sex,
levels = c(1, 2),
labels = c("Hombre", "Mujer")),
ph.ecog2 = factor(ph.ecog,
levels = c(0,1,2),
labels = c("asymptomatic",
"symptomatic but completely ambulatory",
"in bed <50% of the day")))
head(lung)
summary(lung)
inst time status age sex ph.ecog ph.karno
Min. : 1.00 Min. : 5.0 Min. :1.000 Min. :39.00 Hombre:102 Min. :0.0000 Min. : 50.00
1st Qu.: 3.00 1st Qu.: 175.2 1st Qu.:1.000 1st Qu.:57.00 Mujer : 64 1st Qu.:0.0000 1st Qu.: 70.00
Median :11.00 Median : 268.5 Median :2.000 Median :63.50 Median :1.0000 Median : 80.00
Mean :10.69 Mean : 311.1 Mean :1.717 Mean :62.52 Mean :0.9458 Mean : 82.17
3rd Qu.:15.00 3rd Qu.: 422.8 3rd Qu.:2.000 3rd Qu.:70.00 3rd Qu.:1.0000 3rd Qu.: 90.00
Max. :32.00 Max. :1022.0 Max. :2.000 Max. :82.00 Max. :2.0000 Max. :100.00
pat.karno meal.cal wt.loss status2 surv_obj.time surv_obj.status
Min. : 30.00 Min. : 96.0 Min. :-24.000 Min. :0.0000 Min. : 5.0000 Min. :0.0000000
1st Qu.: 70.00 1st Qu.: 616.0 1st Qu.: 0.000 1st Qu.:0.0000 1st Qu.: 175.2500 1st Qu.:0.0000000
Median : 80.00 Median : 975.0 Median : 6.500 Median :1.0000 Median : 268.5000 Median :1.0000000
Mean : 79.64 Mean : 928.2 Mean : 9.657 Mean :0.7169 Mean : 311.0904 Mean :0.7168675
3rd Qu.: 90.00 3rd Qu.:1168.8 3rd Qu.: 15.000 3rd Qu.:1.0000 3rd Qu.: 422.7500 3rd Qu.:1.0000000
Max. :100.00 Max. :2600.0 Max. : 68.000 Max. :1.0000 Max. :1022.0000 Max. :1.0000000
ph.ecog2
asymptomatic :47
symptomatic but completely ambulatory:81
in bed <50% of the day :38
#| echo: true
print(lung$surv_obj)
[1] 455 210 1022+ 310 361 218 166 170 567 613 707 61 301 81 371 520 574 390 12
[20] 473 26 107 53 814 965+ 93 731 460 153 433 583 95 303 519 643 765 53 246
[39] 689 5 687 345 444 223 60 163 65 821+ 428 230 840+ 305 11 226 426 705 363
[58] 176 791 95 196+ 167 806+ 284 641 147 740+ 163 655 88 245 30 477 559+ 450 156
[77] 529+ 429 351 15 181 283 13 212 524 288 363 199 550 54 558 207 92 60 551+
[96] 293 353 267 511+ 457 337 201 404+ 222 62 458+ 353 163 31 229 156 291 179 376+
[115] 384+ 268 292+ 142 413+ 266+ 320 181 285 301+ 348 197 382+ 303+ 296+ 180 145 269+ 300+
[134] 284+ 292+ 332+ 285 259+ 110 286 270 225+ 269 225+ 243+ 276+ 135 79 59 240+ 202+ 235+
[153] 239 252+ 221+ 185+ 222+ 183 211+ 175+ 197+ 203+ 191+ 105+ 174+ 177+
Evaluación inicial: ¿sex viola PH?
mod1 <- coxph(surv_obj ~ sex + age, data = lung)
summary(mod1)
Call:
coxph(formula = surv_obj ~ sex + age, data = lung)
n= 166, number of events= 119
coef exp(coef) se(coef) z Pr(>|z|)
sexMujer -0.43909 0.64462 0.19791 -2.219 0.0265 *
age 0.01679 1.01694 0.01086 1.546 0.1221
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
exp(coef) exp(-coef) lower .95 upper .95
sexMujer 0.6446 1.5513 0.4374 0.9501
age 1.0169 0.9833 0.9955 1.0388
Concordance= 0.597 (se = 0.032 )
Likelihood ratio test= 8.43 on 2 df, p=0.01
Wald test = 8.05 on 2 df, p=0.02
Score (logrank) test = 8.17 on 2 df, p=0.02
Evaluemos el supuesto de Proporcionalidad
zph1 <- cox.zph(mod1)
zph1
chisq df p
sex 0.912 1 0.34
age 0.557 1 0.46
GLOBAL 1.296 2 0.52
plot(zph1)


- Si el p-valor en la prueba
cox.zph() para
sex es significativo o el gráfico tiene tendencia →
violación de PH.
Estimador de Kaplan-Meier por sexo
#| echo: true
fit_km <- survfit(surv_obj ~ sex, data = lung)
ggsurvplot(fit_km,
data = lung,
conf.int = TRUE,
xlab = "DÃas", ylab = "Supervivencia estimada",
title = "Curvas Kaplan-Meier por sexo")


ggsurvplot(fit_km,
data = lung,
fun = "cloglog",
conf.int = TRUE,
xlab = "DÃas", ylab = "Supervivencia estimada",
title = "Curvas Kaplan-Meier por sexo")


Estratificación por sex
#| echo: true
mod_strat <- coxph(surv_obj ~ age + strata(sex), data = lung)
summary(mod_strat)
Call:
coxph(formula = surv_obj ~ age + strata(sex), data = lung)
n= 166, number of events= 119
coef exp(coef) se(coef) z Pr(>|z|)
age 0.01600 1.01613 0.01085 1.475 0.14
exp(coef) exp(-coef) lower .95 upper .95
age 1.016 0.9841 0.9948 1.038
Concordance= 0.558 (se = 0.031 )
Likelihood ratio test= 2.23 on 1 df, p=0.1
Wald test = 2.18 on 1 df, p=0.1
Score (logrank) test = 2.18 on 1 df, p=0.1
zph_strat <- cox.zph(mod_strat)
zph_strat
chisq df p
age 0.503 1 0.48
GLOBAL 0.503 1 0.48
plot(zph_strat)

- La función de riesgo base es diferente para hombres
y mujeres.
- La covariable
age tiene el mismo efecto en ambos
estratos.
Comparación de modelos: sin estratificación vs con
estratificación
#| echo: true
# Modelo sin estratificación
AIC(mod1)
[1] 1001.58
# Modelo con estratificación
AIC(mod_strat)
[1] 854.7006
# Ver ambos resúmenes
summary(mod1)
Call:
coxph(formula = surv_obj ~ sex + age, data = lung)
n= 166, number of events= 119
coef exp(coef) se(coef) z Pr(>|z|)
sexMujer -0.43909 0.64462 0.19791 -2.219 0.0265 *
age 0.01679 1.01694 0.01086 1.546 0.1221
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
exp(coef) exp(-coef) lower .95 upper .95
sexMujer 0.6446 1.5513 0.4374 0.9501
age 1.0169 0.9833 0.9955 1.0388
Concordance= 0.597 (se = 0.032 )
Likelihood ratio test= 8.43 on 2 df, p=0.01
Wald test = 8.05 on 2 df, p=0.02
Score (logrank) test = 8.17 on 2 df, p=0.02
#| echo: true
summary(mod_strat)
Call:
coxph(formula = surv_obj ~ age + strata(sex), data = lung)
n= 166, number of events= 119
coef exp(coef) se(coef) z Pr(>|z|)
age 0.01600 1.01613 0.01085 1.475 0.14
exp(coef) exp(-coef) lower .95 upper .95
age 1.016 0.9841 0.9948 1.038
Concordance= 0.558 (se = 0.031 )
Likelihood ratio test= 2.23 on 1 df, p=0.1
Wald test = 2.18 on 1 df, p=0.1
Score (logrank) test = 2.18 on 1 df, p=0.1
anova(mod1, mod_strat, test = "LRT")
Analysis of Deviance Table
Cox model: response is surv_obj
Model 1: ~ sex + age
Model 2: ~ age + strata(sex)
loglik Chisq Df Pr(>|Chi|)
1 -498.79
2 -426.35 144.88 1 < 2.2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
- Se puede comparar el AIC para ver si mejora el ajuste.
- Nota: el coeficiente de
sex desaparece en el modelo
estratificado.
Comparación visual de curvas ajustadas
#| echo: true
library(survminer)
fit_strat <- survfit(mod_strat)
ggsurvplot(fit_strat, data = lung,
legend.title = "Sexo (estratos)",
xlab = "DÃas", ylab = "Supervivencia ajustada")


¿Qué sigue?
#| echo: true
mod_strat2 <- coxph(surv_obj ~ age + sex + ph.karno + pat.karno + meal.cal +
ph.ecog2 + wt.loss, data = lung)
summary(mod_strat2)
Call:
coxph(formula = surv_obj ~ age + sex + ph.karno + pat.karno +
meal.cal + ph.ecog2 + wt.loss, data = lung)
n= 166, number of events= 119
coef exp(coef) se(coef) z Pr(>|z|)
age 9.968e-03 1.010e+00 1.173e-02 0.849 0.39563
sexMujer -5.489e-01 5.776e-01 2.023e-01 -2.713 0.00667 **
ph.karno 2.213e-02 1.022e+00 1.133e-02 1.954 0.05072 .
pat.karno -1.146e-02 9.886e-01 8.579e-03 -1.336 0.18160
meal.cal 2.708e-05 1.000e+00 2.609e-04 0.104 0.91733
ph.ecog2symptomatic but completely ambulatory 6.474e-01 1.911e+00 2.822e-01 2.294 0.02178 *
ph.ecog2in bed <50% of the day 1.454e+00 4.282e+00 4.621e-01 3.147 0.00165 **
wt.loss -1.430e-02 9.858e-01 7.830e-03 -1.827 0.06775 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
exp(coef) exp(-coef) lower .95 upper .95
age 1.0100 0.9901 0.9871 1.0335
sexMujer 0.5776 1.7314 0.3885 0.8587
ph.karno 1.0224 0.9781 0.9999 1.0453
pat.karno 0.9886 1.0115 0.9721 1.0054
meal.cal 1.0000 1.0000 0.9995 1.0005
ph.ecog2symptomatic but completely ambulatory 1.9106 0.5234 1.0989 3.3217
ph.ecog2in bed <50% of the day 4.2822 0.2335 1.7311 10.5931
wt.loss 0.9858 1.0144 0.9708 1.0010
Concordance= 0.651 (se = 0.029 )
Likelihood ratio test= 26.72 on 8 df, p=8e-04
Wald test = 26.35 on 8 df, p=9e-04
Score (logrank) test = 27.34 on 8 df, p=6e-04
zph_strat2 <- cox.zph(mod_strat2)
zph_strat2
chisq df p
age 0.0541 1 0.816
sex 1.4949 1 0.221
ph.karno 5.5965 1 0.018
pat.karno 3.2368 1 0.072
meal.cal 6.6377 1 0.010
ph.ecog2 3.7785 2 0.151
wt.loss 0.0860 1 0.769
GLOBAL 14.3160 8 0.074
plot(zph_strat2)







modelo_completo <- coxph(surv_obj ~ age + strata(sex) + ph.karno + pat.karno + meal.cal +
ph.ecog2 + wt.loss, data = lung)
birthwt.step <- stepAIC(modelo_completo, trace = TRUE)
Start: AIC=848.42
surv_obj ~ age + strata(sex) + ph.karno + pat.karno + meal.cal +
ph.ecog2 + wt.loss
Df AIC
- meal.cal 1 846.43
- age 1 846.91
- pat.karno 1 848.36
<none> 848.42
- ph.karno 1 849.72
- wt.loss 1 850.36
- ph.ecog2 2 853.97
Step: AIC=846.43
surv_obj ~ age + strata(sex) + ph.karno + pat.karno + ph.ecog2 +
wt.loss
Df AIC
- age 1 844.94
<none> 846.43
- pat.karno 1 846.43
- ph.karno 1 847.74
- wt.loss 1 848.37
- ph.ecog2 2 851.99
Step: AIC=844.94
surv_obj ~ strata(sex) + ph.karno + pat.karno + ph.ecog2 + wt.loss
Df AIC
<none> 844.94
- pat.karno 1 845.00
- ph.karno 1 845.86
- wt.loss 1 847.13
- ph.ecog2 2 850.50
#| echo: true
mod_select <- coxph(surv_obj ~ strata(sex) + pat.karno + ph.karno + ph.ecog2 + wt.loss, data = lung)
summary(mod_select)
Call:
coxph(formula = surv_obj ~ strata(sex) + pat.karno + ph.karno +
ph.ecog2 + wt.loss, data = lung)
n= 166, number of events= 119
coef exp(coef) se(coef) z Pr(>|z|)
pat.karno -0.012249 0.987826 0.008474 -1.445 0.14833
ph.karno 0.018725 1.018901 0.011237 1.666 0.09565 .
ph.ecog2symptomatic but completely ambulatory 0.583929 1.793070 0.282687 2.066 0.03886 *
ph.ecog2in bed <50% of the day 1.419825 4.136396 0.464084 3.059 0.00222 **
wt.loss -0.015643 0.984479 0.007915 -1.976 0.04813 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
exp(coef) exp(-coef) lower .95 upper .95
pat.karno 0.9878 1.0123 0.9716 1.0044
ph.karno 1.0189 0.9814 0.9967 1.0416
ph.ecog2symptomatic but completely ambulatory 1.7931 0.5577 1.0303 3.1205
ph.ecog2in bed <50% of the day 4.1364 0.2418 1.6657 10.2719
wt.loss 0.9845 1.0158 0.9693 0.9999
Concordance= 0.636 (se = 0.03 )
Likelihood ratio test= 19.99 on 5 df, p=0.001
Wald test = 20.75 on 5 df, p=9e-04
Score (logrank) test = 21.54 on 5 df, p=6e-04
zph_select <- cox.zph(mod_select)
zph_select
chisq df p
pat.karno 2.4778 1 0.115
ph.karno 5.1171 1 0.024
ph.ecog2 3.7148 2 0.156
wt.loss 0.0602 1 0.806
GLOBAL 6.5404 5 0.257
plot(zph_select)




LS0tCnRpdGxlOiAiUiBOb3RlYm9vayIKb3V0cHV0OiBodG1sX25vdGVib29rCi0tLQoKYGBge3IgcmVzdWx0cz0naGlkZScsIGVjaG89RkFMU0UsIG1lc3NhZ2U9RkFMU0V9CnBhY2thZ2VzIDwtIGMoInN1cnZpdmFsIiwiYXNhdXIiLCJzdXJ2bWluZXIiLCJkcGx5ciIsIk1BU1MiKQppZiAobGVuZ3RoKHNldGRpZmYocGFja2FnZXMsIHJvd25hbWVzKGluc3RhbGxlZC5wYWNrYWdlcygpKSkpID4gMCkgewogIGluc3RhbGwucGFja2FnZXMoc2V0ZGlmZihwYWNrYWdlcywgcm93bmFtZXMoaW5zdGFsbGVkLnBhY2thZ2VzKCkpKSwgCiAgICAgICAgICAgICAgICAgICByZXBvcyA9ICJodHRwOi8vY3Jhbi5yc3R1ZGlvLmNvbSIpCn0Kc2FwcGx5KHBhY2thZ2VzLCByZXF1aXJlLCBjaGFyYWN0ZXIub25seT1UUlVFKQpgYGAKCiMjIERhdGFzZXQgZWplbXBsbzogYGx1bmdgIChwYXF1ZXRlIGBzdXJ2aXZhbGApIAoKYGBge3J9CiN8IGVjaG86IHRydWUKbGlicmFyeShzdXJ2aXZhbCkKZGF0YShjYW5jZXIsIHBhY2thZ2U9InN1cnZpdmFsIikKdGFibGUobHVuZyRwaC5lY29nKQpsdW5nIDwtIGx1bmcgJT4lIAogIG5hLm9taXQoKSAlPiUKICBtdXRhdGUoc3RhdHVzMiA9IGlmZWxzZShzdGF0dXMgPT0gMiwgMSwgMCkpICU+JQogIGZpbHRlcihwaC5lY29nPDMpICU+JQogIG11dGF0ZShzdXJ2X29iaiA9IFN1cnYodGltZSwgc3RhdHVzMiksCiAgICAgICAgIHNleCA9IGZhY3RvcihzZXgsIAogICAgICAgICAgICAgICAgICAgICAgbGV2ZWxzID0gYygxLCAyKSwgCiAgICAgICAgICAgICAgICAgICAgICBsYWJlbHMgPSBjKCJIb21icmUiLCAiTXVqZXIiKSksCiAgICAgICAgIHBoLmVjb2cyID0gZmFjdG9yKHBoLmVjb2csIAogICAgICAgICAgICAgICAgICAgICAgICAgICBsZXZlbHMgPSBjKDAsMSwyKSwgCiAgICAgICAgICAgICAgICAgICAgICAgICAgIGxhYmVscyA9IGMoImFzeW1wdG9tYXRpYyIsIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICJzeW1wdG9tYXRpYyBidXQgY29tcGxldGVseSBhbWJ1bGF0b3J5IiwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAiaW4gYmVkIDw1MCUgb2YgdGhlIGRheSIpKSkKICAKCmhlYWQobHVuZykKYGBgCgpgYGB7cn0Kc3VtbWFyeShsdW5nKQpgYGAKCmBgYHtyfQojfCBlY2hvOiB0cnVlCnByaW50KGx1bmckc3Vydl9vYmopCmBgYAoKIyMgRXZhbHVhY2nDs24gaW5pY2lhbDogwr9zZXggdmlvbGEgUEg/CgpgYGB7cn0KbW9kMSA8LSBjb3hwaChzdXJ2X29iaiB+IHNleCArIGFnZSwgZGF0YSA9IGx1bmcpCnN1bW1hcnkobW9kMSkKYGBgCgojIyMgRXZhbHVlbW9zIGVsIHN1cHVlc3RvIGRlIFByb3BvcmNpb25hbGlkYWQKCmBgYHtyfQp6cGgxIDwtIGNveC56cGgobW9kMSkKenBoMQpwbG90KHpwaDEpCmBgYAoKLSAgIFNpIGVsIHAtdmFsb3IgZW4gbGEgcHJ1ZWJhIGBjb3guenBoKClgIHBhcmEgYHNleGAgZXMgc2lnbmlmaWNhdGl2byBvIGVsIGdyw6FmaWNvIHRpZW5lIHRlbmRlbmNpYSDihpIgKip2aW9sYWNpw7NuIGRlIFBIKiouCgojIyMgRXN0aW1hZG9yIGRlIEthcGxhbi1NZWllciBwb3Igc2V4bwoKYGBge3J9CiN8IGVjaG86IHRydWUKZml0X2ttIDwtIHN1cnZmaXQoc3Vydl9vYmogfiBzZXgsIGRhdGEgPSBsdW5nKQoKZ2dzdXJ2cGxvdChmaXRfa20sCiAgICAgICAgICAgZGF0YSA9IGx1bmcsCiAgICAgICAgICAgY29uZi5pbnQgPSBUUlVFLAogICAgICAgICAgIHhsYWIgPSAiRMOtYXMiLCB5bGFiID0gIlN1cGVydml2ZW5jaWEgZXN0aW1hZGEiLAogICAgICAgICAgIHRpdGxlID0gIkN1cnZhcyBLYXBsYW4tTWVpZXIgcG9yIHNleG8iKQpgYGAKCmBgYHtyfQpnZ3N1cnZwbG90KGZpdF9rbSwKICAgICAgICAgICBkYXRhID0gbHVuZywKICAgICAgICAgICBmdW4gPSAiY2xvZ2xvZyIsCiAgICAgICAgICAgY29uZi5pbnQgPSBUUlVFLAogICAgICAgICAgIHhsYWIgPSAiRMOtYXMiLCB5bGFiID0gIlN1cGVydml2ZW5jaWEgZXN0aW1hZGEiLAogICAgICAgICAgIHRpdGxlID0gIkN1cnZhcyBLYXBsYW4tTWVpZXIgcG9yIHNleG8iKQpgYGAKCiMjIEVzdHJhdGlmaWNhY2nDs24gcG9yIGBzZXhgCgpgYGB7cn0KI3wgZWNobzogdHJ1ZQptb2Rfc3RyYXQgPC0gY294cGgoc3Vydl9vYmogfiBhZ2UgKyBzdHJhdGEoc2V4KSwgZGF0YSA9IGx1bmcpCnN1bW1hcnkobW9kX3N0cmF0KQpgYGAKCmBgYHtyfQp6cGhfc3RyYXQgPC0gY294LnpwaChtb2Rfc3RyYXQpCnpwaF9zdHJhdApwbG90KHpwaF9zdHJhdCkKYGBgCgotICAgTGEgZnVuY2nDs24gZGUgcmllc2dvIGJhc2UgZXMgZGlmZXJlbnRlIHBhcmEgKipob21icmVzKiogeSAqKm11amVyZXMqKi4KLSAgIExhIGNvdmFyaWFibGUgYGFnZWAgdGllbmUgZWwgbWlzbW8gZWZlY3RvIGVuIGFtYm9zIGVzdHJhdG9zLgoKIyMgQ29tcGFyYWNpw7NuIGRlIG1vZGVsb3M6IHNpbiBlc3RyYXRpZmljYWNpw7NuIHZzIGNvbiBlc3RyYXRpZmljYWNpw7NuCgpgYGB7cn0KI3wgZWNobzogdHJ1ZQojIE1vZGVsbyBzaW4gZXN0cmF0aWZpY2FjacOzbgpBSUMobW9kMSkKCiMgTW9kZWxvIGNvbiBlc3RyYXRpZmljYWNpw7NuCkFJQyhtb2Rfc3RyYXQpCgojIFZlciBhbWJvcyByZXPDum1lbmVzCnN1bW1hcnkobW9kMSkKYGBgCgpgYGB7cn0KI3wgZWNobzogdHJ1ZQpzdW1tYXJ5KG1vZF9zdHJhdCkKYW5vdmEobW9kMSwgbW9kX3N0cmF0LCB0ZXN0ID0gIkxSVCIpCmBgYAoKLSAgIFNlIHB1ZWRlIGNvbXBhcmFyIGVsIEFJQyBwYXJhIHZlciBzaSBtZWpvcmEgZWwgYWp1c3RlLgotICAgTm90YTogZWwgY29lZmljaWVudGUgZGUgYHNleGAgZGVzYXBhcmVjZSBlbiBlbCBtb2RlbG8gZXN0cmF0aWZpY2Fkby4KCiMjIENvbXBhcmFjacOzbiB2aXN1YWwgZGUgY3VydmFzIGFqdXN0YWRhcwoKYGBge3J9CiN8IGVjaG86IHRydWUKbGlicmFyeShzdXJ2bWluZXIpCmZpdF9zdHJhdCA8LSBzdXJ2Zml0KG1vZF9zdHJhdCkKCmdnc3VydnBsb3QoZml0X3N0cmF0LCBkYXRhID0gbHVuZywKICAgICAgICAgICBsZWdlbmQudGl0bGUgPSAiU2V4byAoZXN0cmF0b3MpIiwKICAgICAgICAgICB4bGFiID0gIkTDrWFzIiwgeWxhYiA9ICJTdXBlcnZpdmVuY2lhIGFqdXN0YWRhIikKYGBgCgojIyMgwr9RdcOpIHNpZ3VlPwoKCmBgYHtyfQojfCBlY2hvOiB0cnVlCm1vZF9zdHJhdDIgPC0gY294cGgoc3Vydl9vYmogfiBhZ2UgKyBzZXggKyBwaC5rYXJubyArIHBhdC5rYXJubyArIG1lYWwuY2FsICsgCiAgICAgICAgICAgICAgICAgICAgICBwaC5lY29nMiArIHd0Lmxvc3MsIGRhdGEgPSBsdW5nKQpzdW1tYXJ5KG1vZF9zdHJhdDIpCmBgYAoKCmBgYHtyfQp6cGhfc3RyYXQyIDwtIGNveC56cGgobW9kX3N0cmF0MikKenBoX3N0cmF0MgpwbG90KHpwaF9zdHJhdDIpCmBgYAoKCmBgYHtyfQptb2RlbG9fY29tcGxldG8gPC0gY294cGgoc3Vydl9vYmogfiBhZ2UgKyBzdHJhdGEoc2V4KSArIHBoLmthcm5vICsgcGF0Lmthcm5vICsgbWVhbC5jYWwgKyAKICAgICAgICAgICAgICAgICAgICAgIHBoLmVjb2cyICsgd3QubG9zcywgZGF0YSA9IGx1bmcpCmJpcnRod3Quc3RlcCA8LSBzdGVwQUlDKG1vZGVsb19jb21wbGV0bywgdHJhY2UgPSBUUlVFKQpgYGAKCgoKYGBge3J9CiN8IGVjaG86IHRydWUKbW9kX3NlbGVjdCA8LSBjb3hwaChzdXJ2X29iaiB+IHN0cmF0YShzZXgpICsgcGF0Lmthcm5vICsgcGgua2Fybm8gKyBwaC5lY29nMiArIHd0Lmxvc3MsIGRhdGEgPSBsdW5nKQpzdW1tYXJ5KG1vZF9zZWxlY3QpCmBgYAoKCgpgYGB7cn0KenBoX3NlbGVjdCA8LSBjb3guenBoKG1vZF9zZWxlY3QpCnpwaF9zZWxlY3QKcGxvdCh6cGhfc2VsZWN0KQpgYGAK